An Experimental Investigation of the Scaling of Columnar Joints 
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Columnar jointing is a fracture pattern common in igneous rocks in which cracks self-organize into 
a roughly hexagonal arrangement, leaving behind an ordered colonnade. We report observations of 
columnar jointing in a laboratory analog system, desiccated corn starch slurries. Using measure- 
ments of moisture density, evaporation rates, and fracture advance rates as evidence, we suggest 
an advective-diffusive system is responsible for the rough scaling behavior of columnar joints. This 
theory explains the order of magnitude difference in scales between jointing in lavas and in starches. 
We investigated the scaling of average columnar cross-sectional areas due to the evaporation rate, 
the analog of the cooling rate of igneous columnar joints. We measured column areas in experiments 
where the evaporation rate depended on lamp height and time, in experiments where the evapora- 
tion rate was fixed using feedback methods, and in experiments where gelatin was added to vary the 
rheology of the starch. Our results suggest that the column area at a particular depth is related to 
both the current conditions, and hysteretically to the geometry of the pattern at previous depths. 
We argue that there exists a range of stable column scales allowed for any particular evaporation 
rate. 

PACS numbers: 45.70.Qj, 46.50.-l-a 



INTRODUCTION 

Columnar jointing, the spontaneous fracture of rock 
into regular prismatic columns, was first reported in the 
scientific literature in the 17th century Q. Travelers' re- 
ports of the Giant's Causeway in Antrim, Northern Ire- 
land, excited speculation as to what mechanism could 
produce such a startling, seemingly unnatural display of 
order and symmetry. It is now understood that colum- 
nar joints are common in lava flows found worldwide, in 
a wide variety of materials and deposition conditions. 
An example outcrop is illustrated in figure ^ (a). In 
lava, columns can vary in scale between tiny, cm-sized 
colonnades 0, to massive pillars more than 2 m in di- 
ameter. In spite of several centuries of attention, the 
ordering mechanism is still, like many fracture problems, 
not fully understood. 

The qualitative mechanism of joint formation is easily 
described, but questions of what sets the column scale 
and drives ordering are largely open. The process typ- 
ically begins with a pool of lava, recently formed and 
immobile. The lava can be assumed to have an infinite 
horizontal extent, and cools from both upper and lower 
surfaces. Considering only one of the free surfaces, we 
find that a solidification front moves into the melt as it 
cools. Stresses parallel to the front, caused by thermal 
contraction, lead to a series of fractures perpendicular to 
the front. The fractures advance intermittently, with the 
fracture tips confined to a layer only a few cm thick that 
trails closely behind the solidus. Thus, the pattern can 
be thought of as a dynamic 2D fracture network, that 
leaves a columnar record of itself as it advances in the 



third direction. As the cracks advance into the cooling 
rock, their pattern self-organizes into a network that is 
deceptively close to, but not quite, a regular hexagonal 
lattice 01 ■ Further discussion of the geological process 
can be found in refs. 

The pattern is difficult to study in the field, due to the 
requirement of visualizing large landscapes in 3D. How- 
ever, in the absence of strong observational constraints, 
a number of distinct theories have been developed to ex- 
plain columnar jointing. In general, these studies are 
interested in either scaling (refs. |y, [Dl for example) or 
ordering (eg. 0,13,0) > but rarely produce predictions for 
both. Fortunately, an analog columnar fracture process 
can be observed in desiccated starch slurries, a fact that 
was first systematically studied by MiillerjHj|, although 
it had evidently been noticed much earlier llH Il2l |. In 
starches, the subject of this paper, columns are usually 
between 1-10 mm in diameter (2-3 order of magnitude 
smaller than are typical in lava) , as shown in figure ^ 
(b). This laboratory system offers an unusual opportu- 
nity to do modern, controlled experiments that can test 
the various theories, and perhaps shed light on this clas- 
sic problem in geomorphology ^ 1(X, 1^, 2A, 15^ J,6] . 

The dynamic nature of the columnar crack network 
is one of the features that has stimulated renewed in- 
terest in this old problem. It has been suggested that 
there is a close relationship between a number of dif- 
ferent types of fracture patterns in which the fracture 
network is able to evolve 0, 0, 0| • Thin isotropic me- 
dia undergoing brittle fracture due to shrinkage typically 
display fracture networks which have predominantly 90° 
junctions and a mainly four sided polygonal pattern of 
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FIG. 1: (Color online) Columnar jointing in (a) basalt of the 
Columbia Plateau near Banks Lake (95 cm average diameter), 
and (b) in desiccated corn starch. 

cracks [Hill. Evolved patterns, on the other hand, of- 
ten contain mostly 120° junctions and hexagonal units. 
For example, ice-wedge polygons, which cover much of 
the arctic terrain of the Earth, and also appear to cover 
significant areas on Mars [i^, are seasonal thermal con- 
traction cracks in permafrost. The cracks open during 
the winter and heal during summer as water and/or de- 
tritus fills in the fissures |2l| . The annual cycle allows the 
cracks to move slightly between years, and it has been 
observed that young ice-wedge networks contain more 
four sided polygons than older crack networks, which are 
roughly hexagonal 0|. Septarian concretions, roughly 
spherical bodies which have been fractured and filled in 
with another material, also can display similar hexagonal 
joint patterns 23] . 

Columnar jointing itself also turns out to be a remark- 
ably robust phenomenon. Columnar jointing of some 
kind is found in many igneous rock types, such as basalt 
(including lunar basalt |24| ) , rhyolite, andesite and dacite 
(see reference for more detail). It has been reported 
in other geological materials, such as sandstone and coal, 
and in man-made materials such as smelter slag, glass, 
chemically laced vitrified ice, and certain types of dry- 
ing starches 0, 0, 0, |2^ |2^. Some igneous colum- 
nar formations can be very extensive. An example is 
the Columbia River Basalt Group, which covers half of 



Washington state with columnar basalt hundreds of me- 
ters deep li^. A good understanding of scaling could 
be used to extract quantitative information on cooling 
rate or heat fluxes from such igneous formations through 
simple geometric measurements. 

Furthermore, columnar jointing is also interesting as a 
model problem that tests our understanding of fracture in 
general. It is related to other interesting model problems 
such as thin film fracture (which should be the initiat- 
ing pattern at the free surface), and directional drying 
of thin layers, which produces an analogous 2D dynamic 
ordering of cracks. Both of these are subjects of inten- 
sive study 

III III Enn mill, 

and are important in 
many engineering applications of coatings and thin col- 
loidal films. 

In this paper, we will report our recent experiments 
in desiccating starch. Corn starch is just one of sev- 
eral types of starch which form columns, but we have 
found it to be the most amenable to laboratory experi- 
mentation. We have developed several fully 3D visualiza- 
tion techniques to observe the fracture pattern, as well 
as methods to systematically control the drying process 
0. Using these techniques, we focus in this paper on a 
number of tests of the scaling and ordering mechanism of 
the starch columns. We study scaling within individual 
colonnades, and throughout series of colonnades grown 
under a variety of drying situations, and attempt to in- 
terpret the experimental results in terms of some general 
theoretical considerations. We also present observations 
of the coarsening dynamics underlying the evolution of 
the observed column scale, which suggests some insights 
into the ordering mechanism. 

MATERIALS AND METHODS 

Experiments were performed in mixtures of corn starch 
and water. Our methods are similar to those described by 
Miiller ficj and build on our previous work H] . We desic- 
cated slurries of equal weights of dry corn starch (Canada 
brand) and water, in round, flat-bottomed Pyrex dishes 
(90 mm inner radius), under 250 W halogen heat lamps. 
Traces of bleach were added to sterilize the experiments. 
Samples typically took 1-7 days to dry, and were 3-40 mm 
thick. In some experiments, small amounts of gelatin 
were added to the initial slurry, in order to stiffen the 
mixture. The gelatin (Knox brand) was first dissolved in 
boiling water, mixed with cold water, and then added to 
the slurry, which was allowed to set before the experiment 
began. The initial weights of water and starch were kept 
equal in these experiments, but no bleach was added, as 
it tended to react with the gelatin. 

Fracture in drying corn starch occurs in two distinct 
stages |3l ll,4 Ilt3l|. Initially, a few first-generation (or 
primary) cracks, penetrating the full depth of the layer, 
break the starch layer into large pieces. Subsequent, 
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FIG. 2: (Color online) The total mass of the drying sample, vs. drying time, (a) shows an uncontrolled desiccation, with the 
evaporation rate (b). (c) shows a feedback controlled run in which the final drying rate in (d) is made constant. The three stage 
feedback controlled desiccation shown in (e) and (f) was designed to initiate and develop the columnar jointing in a controllable 
way. Arrows point out the kink in the curves that occurs when the sample ceases to dry homogeneously. In (c)-(f) the dotted 
lines and boxes indicate controlled parameters. 



smaller columnar fractures are then observed to prop- 
agate slowly into the layer from the top drying surface. 
The experiment is complete when the fractures reach the 
base of the slurry, and can be seen to appear in the starch 
on the underside of the Pyrex dish. 

After each sample was dry, a set of columns, chosen 
well away from the edges of the dish, was removed. Sam- 
ples of the colonnades were analyzed in cross-section. 
To avoid edge effects, we focused on the inner half of 
the sample's area, and did not include columns adja- 
cent to first-generation cracks 0, 0. Layers of starch 
were successively and destructively removed by sawing 
away the surface of this sub-sample. Images were ob- 
tained of the surface after each cut, and cross-sectional 
areas were measured using Scion Image analysis soft- 
ware. Typically, 100-200 polygons were measured for 
each layer, and, where shown, error bars represent the 
standard error of the measurements of the mean area. 
This destructive sampling technique produced useful col- 
umn measurements at 2-3 mm depth intervals, but was 
too aggressive to use closer than ^5 mm from the dry- 
ing surface, where the structure is too delicate. When 
higher resolution or near-surface observation was desired, 
we used micro-computed x-ray tomography (Micro-CT) 
to image starch colonnades after desiccation. This tech- 
nique provides full 3D volume filling images, as described 
in Ref. Finally, in experiments where only the final 
cross-sectional area was wanted, we simply photographed 
and analyzed the base of the colonnade after removing it 
from its container. When destructive sampling was used, 
this basal cross-section is also included as the terminal 
data point, to show the consistency of the methodology. 

In order to make sense of the scale of the colonnade, 
it is essential to examine the drying process using con- 
trolled variations of the evaporation rate. We studied 



the effect of several feedback control schemes of increas- 
ing sophistication, shown in Fig. El In all cases where 
evaporation rates are given below, experiments were per- 
formed in 60 mm radius dishes. The mass of the slurry 
was measured once per minute by an electronic bal- 
ance. In the simplest experiments, the drying rate of 
the slurry was varied by manually adjusting the height 
of the heat lamps and the level of ventilation between ex- 
periments. The relative heating intensity was measured 
with a photometer. In these experiments, after an ini- 
tial rapid homogeneous drying phase, we observed that 
the evaporation rate decreased smoothly with time, as 
shown in Fig.[2Ia,b). Hereafter, we refer to these as "un- 
controlled" experiments. To more directly control the 
evaporation rate throughout the experiment, we designed 
a computer-controlled feedback loop that allowed us to 
modulate the duty cycle of the heat lamps, and a small 
overhead fan. We used this technique to mimic the pri- 
mary homogeneous stage of the uncontrolled evaporation 
curve as much as possible, using two types of feedback 
algorithm. In one version, we controlled only the final 
evaporation rate, Q, as shown in Fig. [2Ic,d). For these 
experiments, the controlled evaporation begun slightly 
after columnar jointing was initiated. In order to study 
the initiation of the pattern more carefully, we also im- 
plemented a three stage control algorithm, as shown in 
Fig. [2Ie,f). This control scheme began with a period of 
constant, rapid evaporation, followed by a ramp down 
to a controlled final evaporation rate. The control algo- 
rithm has three free parameters, the initial drying rate 
Qi, the final drying rate Q2, and the slope of the tran- 
sition between them, Q. The initial water fraction for 
all runs was 0.50, by weight. The Q ramp was started 
when the water fraction was reduced to 0.33. Columnar 
jointing typically began late in the ramp phase. In all 
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feedback controlled experiments, the control loop termi- 
nated when the sample was fully dried, or when the duty 
cycle reached 100%, when the algorithm was no longer 
able to maintain the desired evaporation rate. 

RESULTS 

As can be seen in Fig.|2Ia,b), we observed a multi-stage 
evaporation curve under uncontrolled conditions. This 
agrees with the results of Toramaru and Matsumoto , 
but contradicts Miiller, who had suggested that the desic- 
cation curve was smooth . We found that the change 
in the uncontrolled drying curve, marked by the arrow in 
Fig. [2Ia), signals a transition between homogeneous and 
frontal drying. Each has its own style of fracture; the 
former accompanies the primary cracks, while the lat- 
ter accompanies columnar joint formation. In the initial 
stage of the desiccation, the slurry dries uniformly, and 
first-generation cracks open. The plumose structure [l^ 
of these cracks shows that they initiate at the upper dry- 
ing surface, and quickly propagate vertically through the 
entire sample. The resulting disordered fracture networks 
scale like those of ordinary thin films, and have been 
studied as examples thereof 0, |^ . In particular, the 
average spacing between first-generation cracks increases 
roughly linearly with sample depth, and the cracks typi- 
cally meet at 90° junctions. 

When the whole sample reaches a water fraction of 
0.30 ± 0.01, by weight, the slurry begins to dry from the 
top down. This moisture content corresponds to the kink 
in the drying curve indicated by the arrows in figure |21 
Columnar jointing is initiated from the drying surface at 
this time and further drying propagates a fracture front 
into the sample. The front is a clearly delineated bound- 
ary between unfractured starch with a consistency of wet 
clay, and starch which has given way to columnar joint- 
ing. The two layers come apart naturally when a partially 
dried sample is removed from its dish, as all the fracture 
tips are confined to a thin interface, and the thickness of 
the columnar layer can be measured directly with a ruler. 

Fig-Elshows the depth to which columnar jointing had 
progressed in the centre of the sample, in two experi- 
ments. Each point represents columnar depth achieved 
in a separate partially completed desiccation. Fig. |3{a) 
used uncontrolled evaporation with lamps 15 cm above 
the dishes, while Fig. Elb) used a fully controlled evapo- 
ration with Qi — A g/h, Q2 ~ I g/h, and Q = 0.5 g/h^. 
Both experiments had a fully dried depth of 33 mm. In 
the uncontrolled run, the fracture front slows down as 
it enters the sample. We fit the position of the fracture 
front to a power law form with three free parameters, 
the amplitude, exponent and the intercept on the time 
axis. The uncontrolled experiment had an exponent of 
0.60 ± 0.02, which implies that the slowing front pene- 
trates slightly faster than normal diffusion. On the other 
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FIG. 3: Propagation of the fracture front into (a) uncon- 
trolled, and (b) controlled starch desiccation experiments, 
with power law fits. Errors (±1 mm) are approximately the 
size of the data points. 

hand, the corresponding fit to the controlled experiment 
gives an exponent of 0.97 ± 0.06, indicating that con- 
trolling the evaporation rate fixes the average speed of 
the fracture front so that it progresses linearly into the 
sample. 

The columnar fracture front that we observed tracks 
the intrusion of a sharp moisture front. By cutting a 
partially dried sample into thin layers, we measured the 
moisture content as a function of depth by subtracting 
the weights of small samples taken at various depths be- 
fore and after baking them dry. Fig. 0] shows the depth 
dependence of the moisture content of a starch sample 32 
h into an uncontrolled desiccation. A similar sharp dry- 
ing front has also been seen in directional drying of 2D 
suspensions '30^, and has been imaged in drying starch 
slurries using MRI Starch taken from below the 

fracture front closely matches the 0.30 moisture content 
which it achieved before the columnar phase of the desic- 
cation. Above the fracture front, there is a steep decline 
in moisture, starting at a sharp drying front. The front 
is probably maintained by flow driven by capillary pres- 
sure; the "wicking" , of water towards the drying surface, 
rather than by solely by diffusion. These effects are likely 
what lead to the non-diffusive power law behavior seen 
in the fracture advance rate. 

To test that transport did in fact occur through the 
bulk, and was not facilitated by evaporation through the 
cracks, two experiments were set up using C0CI2 and 
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FIG. 4: Fractional moisture content of a partially dried starch 
colonnade. The dashed line indicates the depth of columnar 
fracture. 



CuS04 indicators. 6 g of salt were dissolved in 100 g of 
water, mixed with 100 g of dry starch, and the slurries 
were dried. The location of the salt deposition indicated 
where active evaporation was occurring. In both cases, 
the indicators ended up at the upper drying surface of the 
desiccated sample, mostly within the top 2 mm. C0CI2 
is cobalt blue in its anhydrous state (when dry or warm) , 
and red in its hydrated state (under cool, moist condi- 
tions). In the case of C0CI2, the entire sample of dried 
starch was homogeneously tinted light pink (likely due 
to residual moisture bound to the starch), while the sur- 
face was encrusted with a deep blue deposit. In no case 
was salt left on the walls of the cracks, indicating that 
evaporation did not occur through the cracks. 

These observations imply that vapor is not being trans- 
ported up the cracks to the top surface. This result agrees 
with that of Dufresne et aZ.[2d| who observed the drying 
of nanosphere/ water slurries using CARS microscopy and 
concluded that that evaporation from the drying edge 
dominates over evaporation inside the cracks, where the 
humidity is easily saturated. It remains possible, how- 
ever, that cracks can have some enhancing effect on the 
liquid water transport [l^ . 

We have previously reported our initial investigation 
of scaling behavior in starch colonnades. Some of our 
results are summarized in figure [SJ and more details can 
be found in ref. 0. These results allow us to estimate 
that errors in experimental repeatability between uncon- 
trolled runs is, at most, approximately ±8%. The sta- 
tistical uncertainty in each measurement is smaller than 
this, approximately 4%, and is about the size of the data 
points we present. We found these results to be highly 
repeatable, and not dependent on the absolute depth of 
the sample. For example, we performed two otherwise 
identical experiments on 7.5 and 5.5 cm thick slurries, 
and found that the position of a jump in columnar scale. 
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FIG. 5: Average areas of corn starch columns under uncon- 
trolled drying conditions, vs. depth. The grey diamonds show 
the scales reached at the base of 17 identical slurries dried 
under the same, relatively slow conditions, but with different 
initial thicknesses. Micro-CT x-ray observations on a 25 mm 
deep, more quickly dried slurry (open circles) show a non-zero 
average column area near the drying surface. The filled circles 
show measurements from the destructive sampling of a thick, 
rapidly dried slurry that exhibits a sudden jump in average 
column area. 



as shown in Fig.|Slfor the thicker experiment, occurred at 
the same depth in both cases. In general, for uncontrolled 
dr ying , column scale increases, or coarsens, with depth 
[slllOll We showed that two types of coarsening were pos- 
sible: smooth power-law-likc increases, or sharp, almost 
discontinuous changes in column scale. Using data gath- 
ered by 3D x-ray Micro-CT visualization, we observed a 
non-zero limiting column area near the drying surface. 

Building on our previous studies, we have observed 
identical dishes of starch slurries dried under heat lamps 
at various heights, but without feedback control. Using 
a photometer, we measured the light intensity incident 
on the drying surface, which we assumed to be propor- 
tional to the incident heat flux. Samples 6.5 mm and 
13 mm deep were studied, and the mean cross-sectional 
area of the columns was measured at the base of the 
dried slurry. Fig. |H1 shows the mean column scale as the 
lamp height was varied from 12-80 cm. Faster drying 
generally gives rise to smaller columns. For shallow sam- 
ples, we find an approximately inverse relationship be- 
tween column scale and lamp height, in agreement with 
Toramaru and Matsumoto[l5j, who observed, in simi- 
lar experiments with potato starch, that the final cross- 
sectional area and the initial desiccation rate were in- 
versely related. However, a strict inverse proportionality 
is only found for some of the data in the 6.5 mm thick 
samples. 

We also investigated the effects of adding gelatin to 
the initial slurry. This was done to place some quali- 
tative experimental constraints on ball-and spring mod- 
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FIG. 6: Average column areas at the base of the colonnade 
for 13 mm (open circles) and 6.5 mm (filled circles) deep sam- 
ples. The incident power records photometer readings as the 
heat lamp height was varied. The line represents an inverse 
relationship between drying power and cross-sectional area. 
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FIG. 7: Average cross-sectional area at the base of five 17 
mm deep samples, dried slowly with heat lamp 75 cm above. 
Gelatin concentrations are given relative to the initial mass 
of the slurry. 



els of columnar joint formation. Such "discontinuum" 
fracture models \^ are based on networks of linear 
springs with finite breaking strengths. Fig. [T] shows that 
as gelatin is added, the average final column area in- 
creases, all else being equal. This is what would be ex- 
pected from discontinuum models since the added gelatin 
presumably enhances the Young's modulus and fracture 
toughness of the networks. More detailed measurements 
of the rheology of gelatin-starch mixtures, and more real- 
istic modeling IH^I would be required to pursue this idea 
further. 

The experiments just described show that even the 
simplest level of control has some reproducible effect on 
the final column scale. The column scale was still coars- 



ening in the interior of these samples, however. This 
coarsening is linked to the slowing down of the fracture 
front during uncontrolled desiccation, as was shown in 
Fig. Hi,. 

To study the scaling behavior of corn starch columns 
under more steady-state conditions, we used automated 
feedback to fix the final desiccation rate. In the single- 
control scenario, the overhead lamps were given an ini- 
tial duty cycle of 0.5, and a feedback loop was switched 
on when the desiccation rate dropped to some desired 
value, Q. The feedback loop then fixed the desiccation 
rate through the rest of the experiment, as shown in fig- 
ureEl;c-d). In all constant final column scale was 

eventually chosen, and the coarsening that characterized 
the uncontrolled experiments was halted. Figure |S1 shows 
the results of varying Q from 0.9 to 1.5 g/h. It also in- 
cludes one experiment, where Q = 0.7 g/h, but with an 
initial duty cycle of 0.25. We found the results of these 
experiments to be rather surprising - all four samples 
dried under the same initial conditions, but with differ- 
ent final desiccation rates, were indistinguishable within 
errors. The sample that was initiated differently, but sub- 
ject to the same type of feedback during the controlled 
phase, selected a column area that was remarkably dif- 
ferent. 

This puzzling observation suggested that perhaps long- 
term exposure to heat during the initial phase was slowly 
changing the physical properties of the starch slurries. If 
this was true, the slower initiation of the 0.7 g/h sam- 
ple could have resulted in a sample that had a different 
rheology than the others. We were able to eliminate this 
possibility, however, by preparing six identical dishes of 
slurry, placing them all under one lamp, covering all but 
one with a Pyrex lid, and removing the lid from one addi- 
tional sample every day. The cross-sectional areas at the 
base of all six samples were identical within error. Thus, 
it appears that the subsequent columnar fracture behav- 
ior of the samples is not affected by how long they are 
initially heated without drying. We conjecture that all 
samples that have evolved up to the end of the uniform 
drying stage have identical water content and rheology. 

Experiments using the most sophisticated three stage 
control algorithm confirm that the initiation and early 
development of a colonnade affects its scaling throughout 
the subsequent frontal drying regime. In fully controlled 
experiments, the evaporation rate was fixed throughout 
the entire desiccation, as explained in figure|2e-f). Fig.lHI 
shows the results of changing the control parameter Qi, 
leaving Q2 and Q fixed. Two data sets are shown, for 
which coarsening gradually comes to a halt as a final 
scale is selected. These experiments, which had the same 
final evaporation rate Q2, differ in final scale by over a 
factor of two. The results show that, even given identical 
final evaporation rates Q2, columnar desiccation can se- 
lect average cross-sectional areas that vary considerably, 
depending on Qi a-nd Q- We argue below that this re- 
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FIG. 8: Average cross-sectional area of columns, measured 
through destructive sampling, for single-controlled experi- 
ments with a fixed final desiccation rate Q. The initial duty 
cycle of the desiccating lamps on the 0.7 g/h sample is half 
that of the other samples. Representative errors are shown 
on the 0.7 g/h data. 
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FIG. 9: Average cross-sectional area of columns for two slur- 
ries under fully controlled desiccation rates, measured through 
destructive sampling. The filled circles have Qi = 6 g/h, Q2 
— 1 g/h, and Q = 0.5. The open triangles are identical but 
use Qi = 9 g/h. Representative errors are shown on one data 
set. 



fleets the intrinsically hysteretic nature of the fracture 
process. 

By systematically varying Q2 with the other drying 
parameters remaining fixed, we can explicitly study the 
response of average column area to a family of constant 
final evaporation rates. Results for 600 g of initial slurry, 
Qi = 6 g/h and Q ~ 0.5 are shown in Fig.^| The plot- 
ted cross-sectional areas are averages of 4-5 observations 
made from distinct cuts of destructive sampling, all near 
the bases of the colonnades. The colonnades had a fi- 
nal depth of 33 mm, and average column areas remained 
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columns for experiments with Qi = 6 g/h and Q — 0.5, with 
a varying final evaporation rate Q2- 



constant in approximately the lower half of the colon- 
nades. Although there is a strong dependence of area on 
evaporation rate below Q2 = 1 g/h, there is little to no 
variation in area in the range 1 < Q2 < 1-75 g/h. This 
shows that an established columnar pattern with a given 
scale can propagate stably over a range of evaporation 
rates. 

We have discussed three types of evolution that the av- 
erage column area exhibits: average areas can smoothly 
increase with depth (i.e. coarsen), they can change scale 
catastrophically, or they can arrive at some constant, "se- 
lected" scale. However, this list only describes the behav- 
ior of average column areas. To study the underlying dy- 
namics, we examined the evolution of single columns, us- 
ing 3D x-ray Micro-CT tomograms for both uncontrolled 
and fully controlled desiccations |23.j ■ 

Fig. II II shows the cross-sectional area for a few repre- 
sentative individual columns. The mechanism of coars- 
ening is the termination of joints as they intrude into the 
sample, resulting in the merger of two or three columns 
(3|j. When columns merge, which occurs at the discon- 
tinuities in Fig. ^2 we continued to track the offspring 
column. FigurelTTlfa-b) shows that after a merger event, 
the column rapidly shrinks back to average size, giving 
up area to its neighbors. A column that is spawned from 
a vertex also influences its neighbors, growing at their 
expense. In most cases, we found that such a new col- 
umn merges with one of its neighbors relatively quickly, 
although occasionally one survives and grows to average 
size. An example is shown in Fig. Hir e). The dashed 
arrow indicates the creation of a new column at a vertex 
of the column whose area is plotted. When this happens, 
the observed column shrinks as its new neighbor grows, 
until the two columns merge at the solid arrow. In gen- 
eral, changes in an individual column's area appear to 
involve purely local processes, and affect only the areas 
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FIG. 11: Area evolution for selected individual columns in 
uncontrolled coarsening (a-b) and controlled (c) experiments. 
Column merger events are marked by a dark arrow pointing 
to discontinuities in area. In (c), a dashed arrow points to 
the creation of a new column, which grew from of an adjacent 
vertex. 



of its immediate neighbors. 



DISCUSSION 

The range of experiments discussed above present a 
rather complex picture and require careful interpretation. 
If we restrict ourselves to the question of how the average 
column area that is finally selected scales with final desic- 
cation rate, we find that no single, simple answer emerges 
from the data. Taken together, however, the experiments 
make it possible to draw some important general conclu- 
sions about the mechanisms of scale selection in columnar 
joints in both starch and basalt. 



Advection-diffusion 

The dynamics of columnar jointing in starch can only 
be be understood in conjunction with the dynamics of 
the moisture content within the sample. Given the ob- 
servations shown in Figs. I2I4L we suggest that a nonlinear 
diffusion equation such as 

dtC^V ■[D{c)Vc], (1) 

where the hydraulic diffusivity D depends on the water 
density c, is appropriate to model the water transport 
in the bulk of the starch. We can treat transport by 
the cracks, if any, as part of the boundary conditions 
on this equation. This suggestion is in accord with the 
conclusions of Mizuguchi et al., who performed MRI on 
desiccating samples JJjJ. The porous medium equation 
that can be obtained by from Eqn. ^by setting D(c) oc c 
is known to give rise to propagating self-similar solutions 
with parabolic shape |3J| . Analogously, the shape of the 
sharp front shown in Fig. 0] could be explained by the 
presence of a minimum in D{c) at a moisture content of 
0.3, by mass (a water density of cq = 0.33 g/ml). Such a 
sharp front is expected when capillary forces become the 
dominant driving mechanism for water transport, and it 
indicates a shift from flow-limited to evaporation-limited 
dynamics |3fl| . Dufresne et al. [s^ showed that the cross- 
over length of this type of front scales inversely with the 
evaporation rate, analogous to our parameter Q2. This 
length is characterized by the width of the dropoff in 
moisture at the desiccation front, and cannot be more 
than a few mm in our ^exp eriments. Such fronts are typ- 
ical in drying slurries 28, 

We now consider a starch sample in which the crack 
front and the drying front both move at constant speed v 
in the z direction, as in the fully controlled experiments 
with constant evaporation rate. Moving to the co- moving 
frame introduces an advective term into Eqn. 

dtC^vd,c + V ■[D{c)Vc]. (2) 

Around the minimum in D{c) near the drying front, 
we can make the Taylor expansion 

Dic)^Do + D2{c-co)^ + ... . (3) 

The position of the crack tips follows closely behind the 
desiccation front, so that they will be subject to an ef- 
fective diffusivity close to Dq, but the diffusivity will 
increase rapidly away from the position of the crack 
front. Under these conditions, a constant moisture con- 
tent boundary condition leading the front, and a constant 
flux boundary condition trailing the front are appropri- 
ate. These will quickly lead to a time-independent ver- 
sion of the nonlinear advection diffusion equation, which 
assumes the nondimensional form, 

V^c + Pe' dz'C — nonlinear terms, (4) 
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where z' = zjL^ and the hydrauhc Peclet number, 
Pe' = vL/Dq measures the ratio of advective to diffu- 
sive effects. We take the characteristic length scale L to 
be the effective radius of a column. 

The hydraulic Peclet number offers some interesting 
insights into the relative scale of corn starch columns and 
lava columns. For a typical controlled experiment with 
a fixed evaporation rate of 1 g/h, we found a fracture 
velocity oiv = 0.18±0.01 iim/s, as shown by the slope of 
the data in Fig. Eh. The observed average cross-sectional 
area for these conditions was 2.9±0.2 mm^ using the data 
from Fig. ^] However, as shown in Fig. ^ changes in 
the initiation of the drying process {i.e. in Qi) can cause 
this area to be vary between at least 1.5-4 mm^. We 
therefore choose a value of L in the range 0.7-1.1 mm. 

The diffusion constant Dq is the most difficult property 
to measure, especially since we do not know what the ef- 
fect of the cracks is on the water transport. If we make 
the plausible assumption that it is small (as is certainly 
true for the vapor transport |3nl| ) . then we can estimate 
Dq by fitting a square root function to an uncontrolled 
desiccation curve. Assuming that the position of the frac- 
ture front, tracks a particular concentration c/, and 
that the minima in diffusivity near the desiccation front 
dominates the dynamics, we can estimate the position at 
any time using a solution of the diffusion equation, 

XW = V4AJtcrf-i(c//co). (5) 

which describes diffusion from a half-space with a dry 
dirichlet boundary condition. While not strictly valid in 
this case, it is nonetheless a reasonable first estimate of 
Dq. Given the uncertainties in c/ from the moisture con- 
tent curve, we estimate Dq = (1-5) x 10"^ m^/s. This is 
somewhat lower than the estimate of the diffusivity by 
Miiller, who obtained (0.7-2.7) xlQ-^ mVs 0. How- 
ever, Miiller's value was based on measurements of mass 
changes during only the first five hours of evaporation, 
which probably leads to an overestimate. In fact, it is 
likely that the best value of Dq is near the low end of the 
stated range. Combining these observations, we estimate 
that columnar jointing in corn starch occurs at a value 
of Pe' ^ 0.2, within a factor of 2. 

The advection-diffusion model requires only slight 
modification to apply to the columnar jointing of lavas. 
The thermal formulation of eqn. ^is simply 

V^T + Pc a;,'T = 0, (6) 

where T is the temperature and Pe = vL/k. The thermal 
diffusivity n of basalts does not vary by more than a 
factor of two for temperatures < T < 1500°C IHIli, 
so the nonlinear terms in Eqn. ^ have been neglected. 

In the case of lava, it is known that heat transport by 
the cracks makes an important, indeed dominant, contri- 
bution to the total heat transport. The crack surfaces are 
treated as boundary conditions on Eqn.|3 If the surface 



of the crack is maintained at a constant temperature Tc, 
then the crack heat flux is given in dimensionless terms 
by Newton's law of cooling, 

dy,T = -Bi{T~T,), (7) 

where the Biot number Bi = hL/X. Here, A is the ther- 
mal conductivity and /i is a phenomenological heat trans- 
fer coefficient. The coordinate y' is perpendicular to the 
crack surface. It is believed that cooling in fractured lava 
is greatly enhanced b y th e boiling and reflux of ground- 
water in the cracks j3.l39l| . The high efficiency of this con- 
vective heat transfer process compared to that of thermal 
conduction in lava suggests that the Biot number is large. 
It is important to note that, although the overall cooling 
is driven by the boiling of water, the fracture front in 
lava is located in a region at much higher temperature, 
near 900°C 0|. The fracture tips are therefore likely to 
be in a region where ordinary thermal conductivity dom- 
inates, and up to several meters distant from the cooler 
convective zone. 

If, as we previously assumed, the cracks in the starch 
do not significantly enhance the water transport, then 
the corresponding "hydraulic Biot number" is effectively 
zero. 

We can now estimate the Peclet number Pe relevant 
to the formation of columnar joints in lava. 25 years 
of observations of the gradual solidification of the Ki- 
lauea Iki lava lake, formed in 1959, show that an Eqn. |H| 
accurately describes the isotherms in boreholes and the 
position of the solidification front [s^. These observa- 
tions yield measured values k = 5.0 x 10"'' m^/s, and 
V = 6.7 X 10~* m/s. Observations of columnar jointing 
near the Kilauea caldera of Mauna Loa (although not on 
Kilauea Iki) suggest that columnar jointing at this site is 
typically rather large, with diameters in the range 1.5-3 
m 40]. From these data, we find Pe ~ 0.15, in remark- 
able agreement with the value of Pe' we found for the 
starch experiments. However, both values should only 
be considered accurate to within a factor of 2. 

In summary, we find that the shrinkage fronts in both 
cooling lava and desiccating starch are similar, in dimen- 
sionless terms. Both systems are controlled by advection- 
diffusion processes with similar Peclet numbers. They 
differ, however, in several secondary aspects. The trans- 
port in the bulk of the starch is probably by strongly 
nonlinear diffusion, while the heat transport in the bulk 
of the lava is linear. The role of the cracks in trans- 
port likely also differs. In the starch, it appears most 
likely that the cracks do not significantly enhance water 
transport. In the case of lava, groundwater boiling and 
thermal conduction combine to give a rather complicated 
boundary condition on the temperature and its gradient 
at the cracks. In either case, however, the result of the 
transport is the propagation of a sharp shrinkage front. 
The order of magnitude differences between the scales 
of joints in lavas (~lm columns) and starches ('--^l mm 
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columns) must be directly related to the dynamic similar- 
ity of the shrinkage fronts in the two cases. By imposing 
a constant evaporation rate, we can experimentally probe 
the detailed relationship of the fracture scale to the speed 
of the shrinkage front. 



Scale selection, memory, and hysteresis 

It has often been implicitly assumed that there is a sim- 
ple one-to-one relationship between the scale of columnar 
joints and some externally imposed driving force. For 
joints in lava, it is generally agreed that there is some de- 
pendence of the column cross-sectional area on the cool- 
ing rate 0, 0, IE S ■ The primary geological evidence 
for this relationship comes from studies of the crystal tex- 
ture of the rock in columnar layers separated by regions of 
entablature (smaller scale, irregular jointing). The crys- 
tal texture of the columnar layers is coarser, irnplying 
that it cooled more slowly than the entablature [3, ■ 
Further evidence exists in a rough correlation linking the 
average column area to the widths of fracture features 
known as striae, which are believed to be a proxy for 
cooling rate in basalt In starch experiments 

0,0,13, it is also observed that faster evaporation rates 
generally give smaller columns. While suggestive, these 
qualitative results are insufficient to supply any strong 
constraint, much less to imply a one-to-one relationship, 
between cooling or evaporation rates and average col- 
umn area. Nonetheless, the assumption that there exists 
such a relationship is implicit in much of the modeling of 
columnar jointing phenomena 0, 13,^3 ■ 

This assumption is not universally accepted. Miiller, 
although using the one-to-one assumption in his analysis, 
expresses concern over the lack of experimental proof for 
this principle 0. Saliba and Jagla 4^ have produced 
a finite element calculation of a fracture front intruding 
into a cooling solid, and have made the alternate predic- 
tion that the fracture scale is independent of the driving 
temperature gradient — except possibly during initia- 
tion. We believe that our experiments are the first to 
examine the experimental relationship of columnar joint 
scaling to fixed gradients and their associated fixed frac- 
ture advance rates. 

Our controlled experiments show that when the evap- 
oration rate is explicitly fixed, there are some regimes for 
which there is no strong dependence of the selected aver- 
age cross-sectional area on the final evaporation rate. In 
other regimes, there is such a dependence. Furthermore, 
there are ranges of column scales that are attainable via 
different approaches to the same final evaporation rate. 
And finally, there is the observation of smoothly chang- 
ing evaporation rates leading to sudden jumps in column 
scale 

The experiments of Toramaru and Matsumoto sug- 
gested a simple inverse relationship between average area 



and the evaporation rate This result is particu- 

larly promising since we expect the width of the des- 
iccation front to vary as 1/(52- It is well known that 
uniformly dried, homogeneous thin films have a frac- 
ture spacing that is proportional to the layer thickness 
|l8Ll3lU4lll42l |. However, even though some of our experi- 
ments with varying lamp height yield similar results, they 
cannot be fit as easily to a purely inverse relationship. As 
we showed in , the average areas in these experiments 
would probably continue to evolve if the samples were 
thicker. 

All these observations can be made consistent with 
each other if we assume that the columnar joint pattern 
exhibits a form of dynamic stability, leading to hystere- 
sis. If the pattern has a memory, there may exist a range 
of column scales that are stable for any given evaporation 
rate. The scale observed then depends on the evolution- 
ary history of the pattern. 

By suggesting that the pattern of columnar joints has 
a memory, we are claiming nothing more than that the 
fracture network at any given depth is, in general, simi- 
lar to the fracture network at depths near to it. This is 
suggested by the way that column scale changes. For a 
small increase in average column number of dis- 

crete, local column mergers occur, after which the col- 
umn network slowly relaxes by propagating the change 
through local interactions of columns with their neigh- 
bors. Small, rapid changes over the whole pattern are 
not possible, or at least highly constrained by geome- 
try. Thus, a given pattern scale may be stable over some 
range of conditions. Large, sudden changes of average 
column scale may be possible if a stability threshold is 
passed where, on average, it is more favorable to have all 
columns more than twice as large as they currently are. 
Under these conditions, a large scale change can be en- 
acted by a catastrophic event consisting of many column 
mergers, as we saw with the jump in columnar scale in 
Fig. |31 This has the appearance of a global instability of 
the whole pattern. 

The experiments that we report in Fig. |51 imply that a 
range of column scales could be selected under identical 
final evaporation rates. Alternatively, the data presented 
in Figs . IHI and 1 1 01 show that the same columnar scale can 
be chosen over a range of final drying rates, given the 
same drying initiation. These results suggest a boundary 
between column scales that are stable, and those that 
coarsen, for a given evaporation rate. Different drying 
histories could lead to the selection of different scales 
from within this stable range. 

A superficially appealing argument for a uniquely se- 
lected column scale under given external conditions has 
often been sought in energy considerations. It has been 
asserted that the hexagonal geometry of columnar joints 
arises because such a form optimally releases the strain 
energy per unit area, or otherwise optimally reduces the 
free energy of the system |E ■ this view, the evolu- 
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tion of the columns tends toward a perfectly hexagonal 
lattice, or at least to a static pattern correspondingto a 
local free energy minimum with frozen-in disorder ^ |^ . 
Our results suggest a more dynamical interpretation, in 
which a nonunique scale may be selected within a stable 
band, depending on system history. Such a pattern need 
not be perfectly ordered Q or static [33|. 

A quantitative approach to the energetics of fracture 
patterns is provided by Griffith fracture theory In 
Grifhth's theory, fractures can only advance if the rate of 
increase in surface energy of the opening fracture is less 
than or equal to the rate of decrease in strain energy. The 
difference, if any, between these rates represents dissipa- 
tion; fracture is, in general, a thermodynamically irre- 
versible process. A simplifying assumption often made is 
to assume quasi-static, or steady, slow fracture advance, 
in which case the dissipation is presumed to be small 
and the equality in the Griffith's criterion holds. This is 
a rather restrictive and unrealistic assumption, however. 
In fact, both the rate of strain energy generation and dis- 
sipation are proportional to the velocity of the fracture 
front. Griffith theory can therefore only be used to set 
an upper bound on the fracture density, p expressed as 
the crack area per unit volume, given by 



(8) 



where Ae is the difference in the strain energy density 
between wet and dry samples, and 7^ is the energy re- 
quired to open a unit area of crack surface. Thus, in the 
presence of dissipation, energy considerations alone are 
insufficient to require a one-to-one relationship between 
fracture advance rates and column scale. 

Irwin's modifications to Griffith theory ji^ are some- 
what more useful to the problem of columnar jointing. 
Irwin introduced the concept of fracture toughness, and 
showed how the stress intensification in the vicinity of a 
crack implies that it is much easier to lengthen an existing 
crack than it is to nucleate a new one. Although a quan- 
titative calculation of the stress intensity factor would be 
difficult for a geometry as complex as that of columnar 
joints, the existence of this factor can supply the effect 
of memory to the columnar structure. For a wide range 
of stresses, the downwards extension of existing cracks 
will be favored over the nucleation of new ones, or the 
termination of old ones. Thus, if a network of cracks 
with a given average area is advancing under changing 
conditions, their propagation may be stable, even if the 
changed conditions would not initiate a pattern of that 
scale. The discontinuous jumps we have observed could 
indicate the crossing of a stability boundary between re- 
gions of relatively unchanging crack propagation. 



2D columnar analogs 

Since no comprehensive theory exists for the 3D case of 
columnar joints, it is useful to consider better-studied 2D 
analogs. In these cases, propagating arrays of periodic, 
parallel cracks are observed in behind shrinkage fronts 
moving through thin samples. Both thermal and desic- 
cation experiments have been performed. Crack arrays 
form when a thin layer of wet slurry is dried from one 
edge, both in layers with one free surfac e ll^ls^. a nd in 
suspensions confined between two plates [28, 2^, 3^. The 
thermal experiments typically consist of dipping a heated 
strip of glass, or other material, into a cooler water bath 
at a fixed rate [H^, |^ ■ These 2D fracture processes have 
also been studied theoretically in some detail ,8^^^. 

Boeck et al. have proposed a theory for the case 
of a thermally-quenched strip. Beginning with Eqns. |H| 
and[3for the planar temperature field, and the appropri- 
ate equations for the temperature-dependent stress field, 
they assume quasi-static fracture advance with speed v 
and crack spacing 2L. Under these assumptions, they 
find a one-dimensional manifold of solutions linking L 
and V. Thus, for a given v, there may exist a unique 
spacing L which is stationary in the co-moving frame, 
or there may not. In general, what is actually observed 
depends on stability considerations. Some stationary so- 
lutions will be unstable to small perturbations, leading 
either to new stable solutions with different v and L, or 
possibly to no stable solution. For a restricted set of per- 
turbations, they analyzed a period doubling instability 
where every other crack falls behind the front and pre- 
sumably stops. They showed that there exists a critical 
Peclet number Pcc such that for Pe > Pcc, any value of 
the crack spacing is stable [H^l- depends only very 
weakly on Bi, varying only by a factor of two when Bi 
is varied over four orders of magnitude. While these cal- 
culations cannot be directly applied to the 3D case of 
columnar joints, they do suggest that all columns larger 
than some critical size might be stable for given condi- 
tions, and hence that similar stability criteria might drive 
the observed coarsening towards larger column scales. 
The period doubling instability that Boeck et al. con- 
sidered bears some resemblance to the sudden changes in 
column scale we observed in rapidly dried starch [3| . 

In another calculation, Jagla Q considered the mech- 
anism by which crack arrays in 2D regularize their spac- 
ing under the action of the stress distribution. A pair of 
cracks whose spacing is wider than neighboring pairs de- 
velops an asymmetric stress concentration near the crack 
tips. When these cracks advance, they do so in a di- 
rection that reduces the maximum local stress, and thus 
curve inwards slightly, in accordance with the Cotterell- 
Rice principle |5J|. This tends to produce equidistant 
cracks and hence a regular spacing. The local equaliza- 
tion of column areas that we have observed in Individ- 
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ual columns is presumably due to the action of analo- 
gous stress distributions, although their 3D morphology 
is much more complicated. 

2D experiments also suggest that a range of crack spac- 
ing can be stable for the same external conditions. In 
2D, crack spacin g ge nerally increases linearly with sam- 
ple thickness [lMl36l|. Ronsin and Perrin observed a 
window of allowed crack spacings for a given dipping ve- 
locity in experiments with glass . Similarly, Shorlin et al. 
[isj performed experiments on directionally dried slurries 
on substrates that incorporated a discontinuous step. In 
certain cases, the original pattern of parallel cracks con- 
tinued with the original spacing, even over step which 
doubled the layer thickness. In other cases, the cracks 
responded by doubling their spacing. 

Thus, a variety of experiments and theory in 2D make 
it plausible that a hysteretic, non-unique column scale 
exists for columnar joints, within some band of scales 
limited by stability boundaries. The task then is to locate 
those boundaries and understand the stability properties 
of the evolving fracture networks that lead to columnar 
structures. 



CONCLUSION 

We have studied the scaling behavior of columnar 
joints using three-dimensional observations of desiccated 
starch colonnades. Our results are primarily experimen- 
tal; we measured the evolution of the average column 
area and the size of individual columns under a range 
of controlled and uncontrolled conditions. We found 
a relatively complicated relationship between the desic- 
cation rates and the average column area. Individual 
columns rapidly evolved toward the average area, even 
while coarsening and merging dynamics continued. We 
argued that the regularization of column areas and the 
stress intensity factor in the vicinity of existing cracks 
produces a kind of geometric memory, which allows scale 
selection to be history dependent. We demonstrated this 
most clearly in the fully controlled experiments. In some 
cases there was a selected scale that was stable over a 
range of final evaporation rates. In other cases, where 
the final evaporation rate was fixed, changes in the initi- 
ation influenced the final scale that was selected. These 
results imply that there is a window of column scales that 
are allowed under any given externally applied evapora- 
tion rate. 

We have sketched a theory of the advective concen- 
tration dynamics in which the water moves by nonlinear 
diffusion. By estimating the hydraulic Peclet number for 
this process, and comparing it to the thermal Peclet num- 
ber for columnar jointing in Hawaiian lava lakes, we have 
shown a dynamic similarity that may explain the relative 
scale between columnar jointing in starches and in lavas. 
To complete this theoretical framework, a theory of the 



rheology and fracture properties of the starch needs to 
be formulated, and compared to the corresponding the- 
ory for lava. 

We argued that the scale of columnar joints is proba- 
bly not fully deducible from the energy balance between 
the elastic energy released and the energy consumed to 
create new fracture surfaces. Rather, the column scale 
and its dynamical evolution will require a more complete 
understanding of the stability and instability of propagat- 
ing fracture networks with their associated strain fields 
and their diffusing moisture or temperature fields. Such 
a theory, along with further controlled experiments on 
starch desiccation, as well as new field observations of 
lava columns, should enable some progress toward un- 
ravelling the puzzle of these strangely shaped rocks. 
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